clear;
[bParams, oParams, setts] = model_settings(2);
load([setts.dir,'static_results.mat']);

%% Baseline model solution
tic
[bParams, oParams, setts] = model_settings(2);
[Mom_baseline, dataSim_baseline, f_baseline, entrant, lambda_e, equity_mat, debt_mat, lambda_mat, coupon_mat, optimal_firm, optimal_coupon] = moments(bParams,oParams,setts);
Mom_baseline.tax_benefit = mean(oParams.pi.*dataSim_baseline.leverage(:));
toc
Mom_baseline

%% CS - gamma
tic
[bParams, oParams, setts] = model_settings(2);
bParams.gamma = 1/(1/bParams.gamma*0.95);
[Mom_gamma_low, dataSim_gamma_low, f_gamma_low, ~, ~, equity_mat_gamma_low, ~, ~, ~, ~, optimal_coupon_gamma_low] = moments(bParams,oParams,setts);
Mom_gamma_low.tax_benefit = mean(oParams.pi.*dataSim_gamma_low.leverage(:));
toc

tic
[bParams, oParams, setts] = model_settings(2);
bParams.gamma = 1/(1/bParams.gamma*1.05);
[Mom_gamma_high, dataSim_gamma_high, f_gamma_high, ~, ~, equity_mat_gamma_high, ~, ~, ~, ~, optimal_coupon_gamma_high] = moments(bParams,oParams,setts);
Mom_gamma_high.tax_benefit = mean(oParams.pi.*dataSim_gamma_high.leverage(:));
toc

[Mom_gamma_low; Mom_gamma_high]

%% CS - beta_i 
tic
[bParams, oParams, setts] = model_settings(2);
bParams.beta = bParams.beta*0.8925;
[Mom_betai_low, dataSim_betai_low, f_betai_low, ~, ~, equity_mat_betai_low, ~, ~, ~, ~, optimal_coupon_betai_low] = moments(bParams,oParams,setts);
Mom_betai_low.tax_benefit = mean(oParams.pi.*dataSim_betai_low.leverage(:));
toc

tic
[bParams, oParams, setts] = model_settings(2);
bParams.beta =bParams.beta*1.1075;
[Mom_betai_high, dataSim_betai_high, f_betai_high, ~, ~, equity_mat_betai_high, ~, ~, ~, ~, optimal_coupon_betai_high] = moments(bParams,oParams,setts);
Mom_betai_high.tax_benefit = mean(oParams.pi.*dataSim_betai_high.leverage(:));
toc

[Mom_betai_low; Mom_betai_high]

%% CS - eta
tic
[bParams, oParams, setts] = model_settings(2);
bParams.eta = 0.946;
[Mom_eta_low, dataSim_eta_low, f_eta_low] = moments(bParams,oParams,setts);
Mom_eta_low.tax_benefit = mean(oParams.pi.*dataSim_eta_low.leverage(:));
toc

tic
[bParams, oParams, setts] = model_settings(2);
bParams.eta =1;
[Mom_eta_high, dataSim_eta_high, f_eta_high] = moments(bParams,oParams,setts);
Mom_eta_high.tax_benefit = mean(oParams.pi.*dataSim_eta_high.leverage(:));
toc

[Mom_eta_low; Mom_eta_high]


%% CS - H
tic
[bParams, oParams, setts] = model_settings(2);
bParams.H = bParams.H*0.95;
[Mom_H_low, dataSim_H_low, f_H_low, ~, ~, equity_mat_H_low, ~, ~, ~, ~, optimal_coupon_H_low] = moments(bParams,oParams,setts);
Mom_H_low.tax_benefit = mean(oParams.pi.*dataSim_H_low.leverage(:));
toc

tic
[bParams, oParams, setts] = model_settings(2);
bParams.H = bParams.H*1.05;
[Mom_H_high, dataSim_H_high, f_H_high, ~, ~, equity_mat_H_high, ~, ~, ~, ~, optimal_coupon_H_high] = moments(bParams,oParams,setts);
Mom_H_high.tax_benefit = mean(oParams.pi.*dataSim_H_high.leverage(:));
toc

[Mom_H_low; Mom_H_high]

%% CS - theta
tic
[bParams, oParams, setts] = model_settings(2);
bParams.theta = 0.105;
[Mom_theta_low, dataSim_theta_low, f_theta_low, ~, ~, equity_mat_theta_low, ~, ~, ~, ~, optimal_coupon_theta_low] = moments(bParams,oParams,setts);
Mom_theta_low.tax_benefit = mean(oParams.pi.*dataSim_theta_low.leverage(:));
toc

tic
[bParams, oParams, setts] = model_settings(2);
bParams.theta = 0.145;
[Mom_theta_high, dataSim_theta_high, f_theta_high, ~, ~, equity_mat_theta_high, ~, ~, ~, ~, optimal_coupon_theta_high] = moments(bParams,oParams,setts);
Mom_theta_high.tax_benefit = mean(oParams.pi.*dataSim_theta_high.leverage(:));
toc

[Mom_theta_low; Mom_theta_high]

%% CS - pi
tic
[bParams, oParams, setts] = model_settings(2);
oParams.pi = 0.1;
[Mom_pi_low, dataSim_pi_low, f_pi_low, ~, ~, equity_mat_pi_low] = moments(bParams,oParams,setts);
Mom_pi_low.tax_benefit = mean(oParams.pi.*dataSim_pi_low.leverage(:));
toc

tic
[bParams, oParams, setts] = model_settings(2);
oParams.pi = 0.2;
[Mom_pi_high, dataSim_pi_high, f_pi_high, ~, ~, equity_mat_pi_high, ~, ~, ~, ~, optimal_coupon_pi_high] = moments(bParams,oParams,setts);
Mom_pi_high.tax_benefit = mean(oParams.pi.*dataSim_pi_high.leverage(:));
toc

[Mom_pi_low; Mom_pi_high]

%% CS - xi
tic
[bParams, oParams, setts] = model_settings(2);
bParams.xi = bParams.xi*2;
[Mom_xi_low, dataSim_xi_low, f_xi_low, ~, ~, equity_mat_xi_low, ~, ~, ~, ~, optimal_coupon_xi_low] = moments(bParams,oParams,setts);
Mom_xi_low.tax_benefit = mean(oParams.pi.*dataSim_xi_low.leverage(:));
toc

tic
[bParams, oParams, setts] = model_settings(2);
bParams.xi = bParams.xi*3;
[Mom_xi_high, dataSim_xi_high, f_xi_high, ~, ~, equity_mat_xi_high, ~, ~, ~, ~, optimal_coupon_xi_high] = moments(bParams,oParams,setts);
Mom_xi_high.tax_benefit = mean(oParams.pi.*dataSim_xi_high.leverage(:));
toc

[Mom_baseline; Mom_xi_low; Mom_xi_high]

save([setts.dir,'cs_table.mat']);
load([setts.dir,'cs_table.mat']);

load([setts.dir,'nodebt_results.mat']);
load([setts.dir,'static_results.mat']);

Mom_nd.tax_benefit = mean(oParams.pi.*dataSim_nd.leverage(:));
Mom_ds.tax_benefit = mean(oParams.pi.*dataSim_ds.leverage(:));

%% Build the actual table
CS_tab1 = [Mom_nd; Mom_ds;];
CS_tab2 = [Mom_gamma_high; Mom_baseline; Mom_gamma_low; Mom_betai_low; Mom_baseline; Mom_betai_high; Mom_eta_low; Mom_baseline; Mom_eta_high; Mom_H_low; Mom_baseline; Mom_H_high; Mom_pi_low; Mom_baseline; Mom_pi_high; Mom_baseline; Mom_xi_low; Mom_xi_high; Mom_theta_low; Mom_baseline; Mom_theta_high];

CS_tab1.Properties.RowNames = {'No debt','Static debt'};
CS_tab2.Properties.RowNames = {'$\gamma = 0.475$','$\gamma = 0.500$','$\gamma = 0.525$', ... 
								'$\beta = 3.629$', '$\beta = 4.066$', '$\beta = 4.503$', ...
								'$\eta = 0.946$', '$\eta = 0.973$', '$\eta = 1.00$', ...
								'$H = 3.744$', '$H = 3.941$', '$H = 4.138$', ...
								'$\pi = 0.10$', '$\pi = 0.15$', '$\pi = 0.20$', ...
								'$\phi = 1.09\%$', '$\phi = 2.18\%$', '$\phi = 3.27\%$', ...
								'$\theta = 0.105$', '$\theta = 0.125$', '$\theta = 0.145$'};
								
CS_tab1 = CS_tab1(:,[4 5 6 7 11]);
CS_tab2 = CS_tab2(:,[4 5 6 7 11]);
CS_tab1{:,1:5} = CS_tab1{:,1:5}*100;
CS_tab2{:,1:5} = CS_tab2{:,1:5}*100;

names1 = {'1','No debt','No debt','%.2f';
          '1','Static debt','Static debt','%.2f';};

names2 = {
    '1','$\gamma = 0.475$','$\gamma = 0.475$','%.2f';
    '1','$\gamma = 0.500$','$\gamma = 0.500$','%.2f';
    '1','$\gamma = 0.525$','$\gamma = 0.525$','%.2f';
    '2','','','';
    '1','$\beta = 3.629$','$\beta = 3.629$','%.2f';
    '1','$\beta = 4.066$','$\beta = 4.066$','%.2f';
    '1','$\beta = 4.503$','$\beta = 4.503$','%.2f';
    '2','','','';
    '1','$\eta = 0.946$','$\eta = 0.946$','%.2f';
    '1','$\eta = 0.973$','$\eta = 0.973$','%.2f';
    '1','$\eta = 1.00$','$\eta = 1.00$','%.2f';
    '2','','','';
    '1','$H = 3.744$','$H = 3.744$','%.2f';
    '1','$H = 3.941$','$H = 3.941$','%.2f';
    '1','$H = 4.138$','$H = 4.138$','%.2f';
    '2','','','';
    '1','$\pi = 0.10$','$\pi = 0.10$','%.2f';
    '1','$\pi = 0.15$','$\pi = 0.15$','%.2f';
    '1','$\pi = 0.20$','$\pi = 0.20$','%.2f';
    '2','','','';
    '1','$\phi = 1.09\%$','$\phi = 1.09\%$','%.2f';
    '1','$\phi = 2.18\%$','$\phi = 2.18\%$','%.2f';
    '1','$\phi = 3.27\%$','$\phi = 3.27\%$','%.2f'
    '2','','','';
    '1','$\theta = 0.105$','$\theta = 0.105$','%.2f';
    '1','$\theta = 0.125$','$\theta = 0.125$','%.2f';
    '1','$\theta = 0.145$','$\theta = 0.145$','%.2f'
    };
header = {'','\makecell{Leverage\\Mean}', '\makecell{Leverage\\Var}', '\makecell{Innovation\\Incidence}', '\makecell{Turnover\\Rate} ', '\makecell{Tax\\benefit}'};

ltable(CS_tab1,names1,header,os_dir(setts.dir,'tables\res_table.tex'))
ltable(CS_tab2,names2,header,os_dir(setts.dir,'tables\cs_table2.tex'))